Costurando um perceptron à mão

Para aquecer …

• Vocês leram … 🤭

  • GOODFELLOW, Ian; BENGIO, Yoshua; COURVILLE, Aaron. Deep Learning. Cambridge: MIT Press, 2016. Disponível em: http://www.deeplearningbook.org.

• GitLab

• Todo mundo gosta de quiz..

Modelos

Representação simplificada da realidade

Regressão linear

Mapeia entrada -> saída

Regressão linear

E a regressão logística?

• Já incorpora uma função de ligação.

• Mapeia a entrada para uma probabilidade

• Função de ligação logit \(\phi(x) = \frac{1}{1+e^{-x}}\)

RNA e não linearidade

• Superar as limitações de funções lineares

• Também usa multiplicação de matrizes + funções de ativação não lineares.

  • (Multilayer) perceptron
  • Feedforward

• Aproximar uma função \(y = f(x)\) com \(y = f(x, \theta)\) aprendendo \(\theta\) que melhor aproximação da função

• Requer escolher a forma de otimização, função de custo, e função de saída

Teorema Aproximação Universal

Uma rede neural com uma única camada oculta contendo neurônios suficientes pode aproximar qualquer função contínua em um intervalo fechado, com uma função de ativação adequada.

Representation Learning

• Rede composta de funções \(f(x) =f^{(3)}(f^{(2)}(f^{(1)}(x)))\).

O que vamos representar hoje…

set.seed(1.2025)
m.obs <- 100000
dados <- tibble(x1.obs=runif(m.obs, -3, 3), 
                x2.obs=runif(m.obs, -3, 3)) %>%
  mutate(mu =abs(x1.obs^3 - 30*sin(x2.obs) + 10), 
         y=rnorm(m.obs, mean=mu, sd=1))

Superfície

Arquitetura da rede

Matematicamente

onde \(\phi(x) = \frac{1}{1+e^{-x}}\) representa a função de ativação logística (sigmoide).

\(f(x_{1i}, x_{2i}; \mathbf{\theta})=\hat{y}_i=\phi(x_{1i} w_1 + x_{2i} w_3 + b_1) w_5 + \phi(x_{1i} w_2 + x_{2i} w_4 + b_2) w_6 + b_3.\)

Matricialmente

onde

\(\mathbf{W}^{(1)}=\)

\(\mathbf{W}^{(2)}=\)

\(\mathbf{b}^{(1)}=\)

\(\mathbf{x}=\)

\(\mathbf{h}=\)

\(\mathbf{a}=\)

Feed forward

Feed forward

## sigmoide
phi <- function(x){
  1/(1 + exp (-x))
}

foward_prop <- function(theta, x){
  
W1 <- matrix(c(theta[1], theta[2], theta[3], theta[4]), 2, 2, byrow=T) 

W2 <- matrix(c(theta[5],theta[6]), 2, 1, byrow=T)

b <- c(theta[7], theta[8])

a <- t(W1)%*%x + b 
h <- phi(a)
y_hat <-t(W2)%*%h + theta[9]
return(as.numeric(y_hat)) }

theta_inicial <- rep(0.1,9)
obs <- c(1,1) # o vetor assim já é tratado como coluna

foward_prop(theta_inicial,obs)

[1] 0.2148885

Função de perda

\[ J(\mathbf{\theta}) = \frac{1}{m} \sum_{i=1}^m L(f(x_{1i}, x_{2i}; \mathbf{\theta}), y_i) = \frac{1}{m} \sum_{i=1}^m (y_i - \hat{y}_i)^2, \] onde \(x_{ji}\) representa a \(j\)-ésima covariável (feature) da \(i\)-ésima observação, \(\mathbf{\theta} = (w_1, \ldots, w_6, b_1, b_2, b_3)\) é o vetor de pesos (parâmetros).

custo <- function(y, y_hat){
  mean((y-y_hat)^2)
}

Função de perda

Split

n_treino <- 0.8*m.obs # proporcao da base que é treino

covariaveis_treino <- as.matrix(dados[1:n_treino,1:2])
covariaveis_teste <- as.matrix(dados[(n_treino+1):m.obs,1:2])

observacao_treino <- pull(dados[,4])[1:n_treino]
observacao_teste<- pull(dados[,4])[(n_treino+1):m.obs]

previsoes_teste <- foward_prop(theta_inicial, 
                            t(covariaveis_teste)) # tem q transpor 


custo1 <- custo(observacao_teste, previsoes_teste)
custo1

[1] 659.6568

Gradiente descendente

• Não tem convergência garantida quando aplicada a funções não convexas.

• Sensível aos valores de inicialização

• Nas redes, estamos sempre usando o descendente de uma forma ou de outra …

• Treinar pelo MSE para distribuições Gaussianas, equivale à máxima verossimilhança

Gradiente descendente

Mínimos locais e globais

Inicialização

Momentum

Back propagation

• Durante o treinamento, a propagação direta continua até produzir um escalar de custo ( \(J(\theta)\) ).
  

• Esse custo representa o erro da rede ao comparar a previsão com o valor real.

• O algoritmo de retropropagação (backprop, Rumelhart et al., 1986) permite que a informação do custo flua para trás pela rede.
  

• O objetivo é calcular o gradiente da função de custo em relação aos parâmetros da rede.
  

Back propagation

• A equação analítica do gradiente é simples, mas sua avaliação numérica pode ser computacionalmente cara.
  

• O backpropagation oferece um procedimento simples e eficiente para esse cálculo.
  

• Nos algoritmos de aprendizado, o gradiente mais relevante é \(\nabla_{\theta} J(\theta)\), ou seja, a derivada da função de custo em relação aos parâmetros da rede.

Back propagation: mitos

• Não é um algoritmo de aprendizado completo: Ele apenas calcula o gradiente, enquanto outro algoritmo (ex.: descida do gradiente estocástica) usa esse gradiente para atualizar os pesos.
  

• Não é exclusivo de redes neurais profundas: Ele pode ser aplicado para calcular derivadas de qualquer função, desde que a derivada seja bem definida.
  

Regra da cadeia

• A regra da cadeia no cálculo é usada para calcular derivadas de funções compostas.
  • O backpropagation a aplica para calcular o gradientes.
• Segue uma ordem específica de operações para garantir eficiência computacional e evitar a explosão exponencial de cálculos redundantes.
• O algoritmo de backpropagation é projetado para reduzir o número de subexpressões repetidas, sem priorizar o uso de memória.

Back - propagation

Use a regra da cadeia para encontrar expressões algébricas para o vetor gradiente.\[\nabla_\theta J(\theta) = \left(\frac{\partial J}{\partial w_1}, \ldots, \frac{\partial J}{\partial b_3} \right)\]

\[\nabla J_\theta = \sum_{i=1}^m(\frac{\partial J}{\partial \hat{y}_i} \nabla_\theta \hat{y}_i)\]

Então, para cada observação, temos que o primeiro elemento desse produto é o escalar:

\[\frac{\partial J}{\partial \hat{y}_i} = -2 ( y_i - \hat{y_i})\] Já o segundo elemento desse produto é, no caso dessa rede, um vetor de tamanho 9 onde:

Para o nono elemento (terceiro viés):

\[\frac{\partial \hat{y_i}}{\partial b_3} = 1\]

Para os pesos da cama de saída: \[\frac{\partial \hat{y_i}}{\partial w_6} = \frac{1}{1+e^{-a_2}}\] \[\frac{\partial \hat{y_i}}{\partial w_5} = \frac{1}{1+e^{-a_1}}\]

Onde \(a_1 = x_1w_1 + x_2w_3 + b_1\) e \(a_2 = x_1w_2 + x_2w_4 + b_2\)

Para os vieses da camada intermediária

\[ \frac{\partial \hat{y_i}}{\partial b_2} = w_6\frac{e^{-a_2}}{(1+e^{-a_2})^2} \]

\[ \frac{\partial \hat{y_i}}{\partial b_1} = w_5\frac{e^{-a_1}}{(1+e^{-a_1})^2} \]

Para os pesos entre a camada de entrada e a intermediária:

\[\frac{\partial \hat{y_i}}{\partial w_4} = w_6\frac{e^{-a_2}}{(1+e^{-a_2})^2}x_{2i}\] \[\frac{\partial \hat{y_i}}{\partial w_3} = w_5\frac{e^{-a_1}}{(1+e^{-a_1})^2}x_{2i}\] \[\frac{\partial \hat{y_i}}{\partial w_2} = w_6\frac{e^{-a_2}}{(1+e^{-a_2})^2}x_{1i}\] \[\frac{\partial \hat{y_i}}{\partial w_1} = w_5\frac{e^{-a_1}}{(1+e^{-a_1})^2}x_{1i}\]

Regra da Cadeia no Back propagation

Algoritmo

No R

phi_linha <- function(x){
  exp(-x)/(1+exp(-x))^2}
  

gradiente <- function (theta, x, y ){
  
  n <- length(y)  
  W1 <- matrix(c(theta[1:4]), 2, 2, byrow=T) 
  W2 <- matrix(c(theta[5:6]), 2, 1, byrow=T)
  b <- c(theta[7], theta[8])

  a <- t(W1)%*%t(x) + b 
  h <- phi(a)
  y_hat <- t(W2)%*%h + theta[9]
  y_hat <- as.numeric(y_hat)

  # derivada da perda
  d_Jota<- -2*(y-y_hat)

  # vies 2 camada (saida)
  grad_b3 <- mean(d_Jota *1)

  # pesos 2 camada
  grad_W2camada <- d_Jota * t(h)
  grad_w6 <-  mean(grad_W2camada[,2] )
  grad_w5<- mean(grad_W2camada[,1])

  # primeira camada
  hlinha <- t(phi_linha(a))
  W_2camada_vetorial <- cbind(rep(W2[1], n ),rep(W2[2], n ) )

  # vieses 1 camada
  grad_b <- d_Jota * W_2camada_vetorial * hlinha
  grad_b2 <- mean(grad_b[,2])
  grad_b1 <- mean(grad_b[,1])

  # pesos primeira camada
  grad_W4_3 <- grad_b * x[,2]
  grad_W2_1 <- grad_b * x[,1]

  grad_w4 <- mean(grad_W4_3[,2])
  grad_w3 <- mean(grad_W4_3[,1])
  grad_w2 <- mean(grad_W2_1[,2])
  grad_w1 <- mean(grad_W2_1[,1])


  c(grad_w1, grad_w2, grad_w3,
    grad_w4, grad_w5, grad_w6,
    grad_b1, grad_b2, grad_b3 )
}

Cáculo

theta <- rep(0.1,9)

gradiente1 <- gradiente(theta, covariaveis_treino, observacao_treino)

gradiente1

[1]  -0.1767887  -0.1767887   0.6458047   0.6458047 -22.3246918 -22.3246918
[7]  -1.0732662  -1.0732662 -43.4113383

Fitting

otimizacao <- function(iteracoes=100, # default assim, mas pode mudar
                       theta_inicial=rep(0,9), 
                       x_treino, y_treino, 
                       x_teste=x_treino, 
                       y_teste=y_treino,
                       LR = 0.1){
  
  # não sei se precisava disso, mas achei melhor ir salvando os resultados
  # intermediários em uma lista 
  Perdas_treino<- numeric()
  Perdas_teste <- numeric()
  Perdas <- list()
  Gradientes <- list()
  Previsoes <- list()
  Parametros <- list()
  Parametros[[1]] <- theta_inicial


for (i in 2:(iteracoes+1)){

Gradientes[[i-1]] <- gradiente(Parametros[[i-1]], x_treino , y_treino)
Parametros[[i]] <- Parametros[[i-1]] - LR*(Gradientes[[i-1]])
previsoes_treino <-foward_prop(Parametros[[i-1]], t(x_treino))
previsoes_teste <- foward_prop(Parametros[[i-1]], t(x_teste))
Perdas_treino[i-1] <- custo(y_treino, previsoes_treino)
Perdas_teste[i-1] <- custo(y_teste, previsoes_teste)
#i <- i+1
}
Perdas <- tibble(Perdas_treino,Perdas_teste)
list(Gradientes=Gradientes,theta=Parametros,
     J_treino=Perdas_treino, J_teste=Perdas_teste,
     J= Perdas)
}


## ESSA PARTE É SE QUISER FAZER SÓ EM UM DOS GRUPOS

so_treino <- otimizacao(x_treino=covariaveis_treino,
                        y_treino=observacao_treino )

so_teste <-  otimizacao(x_treino=covariaveis_teste,
                        y_treino=observacao_teste)



# GRADIENTE CALCULADO NO TREINO, MAS PERDA NO TESTE TAMBÉM


treino_teste <- otimizacao(x_treino=covariaveis_treino,
                           y_treino=observacao_treino, 
                           x_teste=covariaveis_teste, 
                           y_teste=observacao_teste  )
minimo <- which(treino_teste$J$Perdas_teste==min(treino_teste$J$Perdas_teste))
theta_ajustado <- treino_teste$theta[[minimo]]

Ao final …

Previsões

Gradiente descendente estocástico

Gradiente descendente estocástico

rapidez <- microbenchmark(
 gradiente(theta, matrix(covariaveis_treino[1:300,],300,2 ), observacao_treino[1:300]),
 gradiente(theta,matrix(covariaveis_treino[1:n_treino,],n_treino,2 ), observacao_treino[1:n_treino])
  )


rapidez

Unit: microseconds
                                                                                                        expr
                gradiente(theta, matrix(covariaveis_treino[1:300, ], 300, 2),      observacao_treino[1:300])
 gradiente(theta, matrix(covariaveis_treino[1:n_treino, ], n_treino,      2), observacao_treino[1:n_treino])
     min       lq      mean   median        uq      max neval cld
   174.3   267.65   519.952   439.15    650.85   5129.7   100  a 
 35749.4 44794.55 89797.907 71320.80 114602.70 384682.1   100   b

Regressão Linear

$Y_i = N(0 + {1}x_{1i}+ 3x{2i}, ) $

• incluíndo mais covariáveis

\(E(Y|x_1, x_2) = \beta_0 + \beta_1x_1 + \beta_2x_2+\beta_3x_1^2 + \beta_4x_2^2 + \beta_5x_1x_2\)

# vendo quem está dentro ou fora
cores_LM <- predict(modelo_linear1,newdata = base_teste,
                 interval="prediction") %>%
  as_tibble() %>% 
  mutate(y_teste=observacao_teste,
         contem = y_teste>=lwr &y_teste<=upr) %>% 
  pull(contem)

base_previsao_teste_RN %>% 
  mutate(ic=cores_LM) %>% 
  ggplot(aes(x1, x2, color=ic))+
  geom_point( size=0.8, alpha=0.7)+
  scale_color_manual("Está no IC", values=c("#FCFC30", "#465EFF"), labels=c("Sim", "Não"))+
  theme_classic()+
  ggtitle(glue("Cobertura = {round(mean(cores_LM),2)}%\n MSE = {round(summary(modelo_linear1)$sigma^2, 2)}"))

cores_RN <- base_previsao_teste_RN %>% 
  mutate(normal=(y-y_hat)/sqrt(J),
         contem = normal>=-1.96 & normal <=1.96) %>% 
  pull(contem)
  
base_previsao_teste_RN %>% 
  mutate(ic=cores_RN) %>% 
  ggplot(aes(x1, x2, color=ic))+
  geom_point( size=0.8, alpha=0.7)+
  scale_color_manual("Está no IC", values=c("#FCFC30", "#465EFF"), labels=c("Sim", "Não"))+
  theme_classic()+
  ggtitle(glue("Cobertura = {round(mean(cores_RN),2)}%\n MSE = {J}"))

Obrigada!